#Part 1 The purpose of this notebook is to determine the selective benefit of glucosinolate and flavonoid compounds through plant competition by creating selection gradients. Selection gradients will involve final body mass as the proxy for fitness, but this will be replaced with fitness once the measurement is in. Concentration will be on the x-axis. This analysis will be done in each treatment seperately and account for family and greenhouse location.
#Part 2 The second purpose is to determine if glucosinolates and flavonoids influence suscpetibility to pathogens and if this susceptibility results in increased fitness
Read in and prep data
library(ggplot2)
library(dplyr)
library(lme4)
#read in data
rm(list=ls())
dat<-read.csv("DataSynthesis.csv")
#Remove maple data from the data set for the time being.
dat<-dat[!grepl("aple",dat$Tag, fixed =T),]
#(including maple controls):
dat<-dat[!dat$treatment=="mcnt",]
#Removing NA's
dat<-dat[!is.na(dat$Sample),]
dat<-dat[!is.na(dat$gluc_Conc),]
#Assign family column
prefamily<-gsub("*.\\|","",dat$Tag)
dat$Family<-gsub("\\-.*","",prefamily)
#remove those with fertilizer treatment, and extra genotypes that are only in the alone treatment.
dat<-dat[!grepl("i",dat$Sample,fixed=T),]
#Correlation between Glucosinolate and Flavonol Concentration
#Average Duplicates.
dat2<-dat2 %>% select(-X) %>% group_by(Tag) %>% summarize(ChlorA=mean(ChlorA),ChlorB=mean(ChlorB),gluc_Conc=mean(gluc_Conc),flav_Conc=mean(flav_Conc),Family=first(Family),GA3=first(GA3),treatment=first(treatment),gh_row=first(gh_row),gh_bench=first(gh_bench),GM_TotalLeaf_Area=first(GM_TotalLeaf_Area),comp_number=first(comp_number),ThripsDam=mean(ThripsDam),WhiteFungDam=mean(WhiteFungDam),BlackPathDam=mean(BlackPathDam),Fern=mean(Fern),gh_col=first(gh_col))
package ‘bindrcpp’ was built under R version 3.4.4
#Standardizing greenhouse location The purpose of this is to get one numeric vector which can be used to determine the location of individuals in the greenhouse. It will summarize the distance to the walls of the greenhouse.
plot(dat2$gluc_Conc~dat2$gh_bench*dat2$gh_row)
#Investigating genetic differences
boxplot(gluc_Conc~Family,data=dat2[dat2$treatment=="a",])
boxplot(gluc_Conc~Family,data=dat2[dat2$treatment=="m",])
boxplot(gluc_Conc~Family,data=dat2[dat2$treatment=="gm",])
boxplot(GM_TotalLeaf_Area~Family,data=dat2[dat2$treatment=="a",])
boxplot(GM_TotalLeaf_Area~Family,data=dat2[dat2$treatment=="m",])
boxplot(GM_TotalLeaf_Area~Family,data=dat2[dat2$treatment=="gm",])
#I will nest as bench/row/collumn and determine from there what effects need to be retained or not via anova.
dat2$GM_TotalLeaf_Area
[1] 8044 9224 8622 12689 10807 3755 12641 6650 12189 5368 9140
[12] 6001 3598 10739 15374 9329 5632 10014 5679 6867 4419 2478
[23] 8360 4515 8532 7050 8540 11180 4570 8811 4884 3500 7054
[34] 8992 9971 9950 NA 6598 10305 8444 3690 10584 11211 5415
[45] 14121 3967 3456 7710 10237 5094 4051 7518 7405 6561 6461
[56] 8124 2024 11096 2912 2988 13042 8765 9731 8220 3558 6068
[67] 6162 12912 4050 11320 4234 12786 4682 5457 3883 6210 11662
[78] 12334 2560 10012 9520 8712 3036 5357 7115 6275 5619 4092
[89] 7113 4083 1260 3065 2365 7518 16025 5094 12801 8684 10321
[100] 11129 8126 13018 5974 11511 13782 11310 2115 4445 10892 4760
[111] 2067 10861 13223 5121 2512 3665 6481 8329 4608 10404 10144
[122] 1736 3480 10920 6780 5988 5647 3010 4075 5648 6410 11052
[133] 2496 8221 9395 4563 8731 4122 3034 9917 10708 6818 10769
[144] 7671 9315 12836 14115 11327 7871 12862 4874 12576 1970 10122
[155] 11602 9984 11558 10532 9761 2304 6574 4122 6356 3725 6434
[166] 7360 4813 8069 9681 9477 7127 10118 4906 3753 6358 7164
[177] 11068 4830 10546 10942 14747 1178 6461 9896 10000 11566 10757
[188] 10358 5601 5934 6281 2520 15967 11153 3862 3092 4796 8456
[199] 9524 9131 9111 12444 8617 5816 3971 5492 5000 4543 7005
[210] 5641 11188 4732 1800 9803 4549 13342 8548 8320 11793 8308
[221] 8092 4920 9979 14448 5970 3477 2646 1592 7780 3448 8475
[232] 12764 9444 2492 10512 4963 12668 11241 7252 11092 9335 6462
[243] 10907 9341 6700 8025 5926 6316 16542 6942 1452 12441 13308
[254] 3612 9806 6969 11863 8976 6572 8411 9304 4889 1748 4812
[265] 7035 13031 4467 13313 2500 3096 6194 972 2122 6370 9504
[276] 4540 6385 10521 7806 15152 8640 10983 9665 8217 4808 10717
[287] 5849 8445 11872 11024 7857 9671 9456 2304 5696 3488 10572
[298] 5992 5434 9953 6930 3854 2550 1691 7227 NA 10215 14991
[309] 4745 3658 7851 14197 9300 4544 10927 4403 8010 8046 8803
[320] 8923 9874 9572 8305 1841 13313 2377 4794 10299 6984 8240
[331] 8309 12474 8522 7805 4149 9960 7687 11250 9343 10524 4928
[342] 11170 5209 7954 6435 7426 13498 7053 8718 6173 5439 4900
[353] 10695 11803 4950 11456 13132 9718 5993 12159 12111 5502 9821
[364] 6174 17238 1923 4768 2088 6471 11503 10236 10909 14955 9739
[375] 12650 2838 5798 5499 7320 16392 7070 NA NA 9558 9050
[386] 8730 9496 12535 9540 10788 11344 11489 7651 6553 1295 7059
[397] 9187 9993 10512 7858 17606 6060 5167 10399 10159 5003 9306
[408] NA 7189 10248 6906 6397 11965 8315 12173 8018 11156 7028
[419] 3306 10267 1800 17244 9232 5616 7809 10319 13392 12559 10927
[430] 918 5112 6342 8499 5472
#Effect of glucosinolates & flavonols on fitness (bodymass for proxy)
#removing those without fitness measurements
dat2<-dat2[!is.na(dat2$GM_TotalLeaf_Area),]
#in the maple treatment
ggplot(dat2[dat2$treatment=="m",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc))
#in the garlic mustard treatment
ggplot(dat2[dat2$treatment=="gm",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc))
#alone
ggplot(dat2[dat2$treatment=="a",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc))
#there does not appear to be any on competition what so ever, however, maybe trends will appear after accounting for gh location. This looks at whether there is an interaction between gluc conc and treatment, which is needed to infer a benefit of glucosinolate concentration.
library(lattice)
PlotResult<-groupedData(GM_TotalLeaf_Area~gluc_Conc|Family,data=dat2[dat2$treatment=="m",])
PlotResult<-groupedData(GM_TotalLeaf_Area~gluc_Conc|Family,data=dat2[dat2$treatment=="gm",])
plot(PlotResult)
#Event though there does not appear to be a relationship for any given genotype, when we account for their average differences across treatment and greenhouse bench, a trend comes out.
#Visualizing mean result.
dat3<-dat2 %>% group_by(Family,treatment) %>% summarize_if(is.numeric,mean)
#in the maple treatment
ggplot(dat3[dat3$treatment=="m",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc,colour=Family))
#in the garlic mustard treatment
ggplot(dat3[dat3$treatment=="gm",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc))
#alone
ggplot(dat3[dat3$treatment=="a",])+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc))
#Remove genotypes that had thier competitor die while in the GM treatment. This needs to be done!!!
table(dat2$comp_number)
1 2 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
1 2 1 2 2 1 1 2 1 1 1 2 2 2 2 1 2 1 2 2 1 1 2 2
26 27 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 51
2 2 1 2 1 2 2 2 2 1 1 2 2 2 1 2 1 1 2 1 1 2 2 2
52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 71 72 73 74 75 76
2 2 2 2 1 2 2 1 1 1 1 2 1 1 2 2 2 2 2 2 1 2 2 2
78 79 80 81 82 83 84 85 86 87 88 89 90 91 92
1 1 2 2 2 2 2 2 2 1 1 1 1 2 2
dat
dat2
#Based on these results, it seems as though i need to include family to avoid psuedo replication, and i need to include bench, but adding greenhouse row does not seem to increase the predictive power that much so i will leave that out. Now i will compare all variables which i am interested in predicting fitness. with these random effects.
#This is the full model, i expect that the influence of gluc conc and flav conc could vary between treatments because of allelopathy, but i have no reason to think that pathogens could influence fitness differently in different treatments. The same goes with ferns. Competition from ferns would affect fitness equally in all treatments, even though there may be more ferns in certain treatments, but this is not what i am testing wtih this analyis. Let us continue with backwards model selection.
#First lets ensure we are using the correct random effects.
fitfull0<-lmer(GM_TotalLeaf_Area~treatment*gluc_Conc*flav_Conc+BlackPathDam+WhiteFungDam+ThripsDam+Fern+(1|Family), data=dat2)
fitfull<-lmer(GM_TotalLeaf_Area~treatment*gluc_Conc*flav_Conc+BlackPathDam+WhiteFungDam+ThripsDam+Fern+(1|Family)+(1|gh_bench), data=dat2)
fitfull2<-lmer(GM_TotalLeaf_Area~treatment*gluc_Conc*flav_Conc+BlackPathDam+WhiteFungDam+ThripsDam+Fern+(1|Family)+(1|gh_bench/gh_row), data=dat2)
fitfull3<-lmer(GM_TotalLeaf_Area~treatment*gluc_Conc*flav_Conc+BlackPathDam+WhiteFungDam+ThripsDam+Fern+(1|Family)+(1|gh_bench/gh_col), data=dat2)
#Gh_Col is the best predictor for the data by far, so i will use that as the random effects in the model.
anova(fitfull0,fitfull) #Evidence to use bench at least.
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
object: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
object: WhiteFungDam + ThripsDam + Fern + (1 | Family)
..1: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
..1: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
object 18 7640.8 7712.9 -3802.4 7604.8
..1 19 7601.0 7677.2 -3781.5 7563.0 41.724 1 1.051e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
anova(fitfull,fitfull2) #No evidence to use row.
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
object: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
object: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench)
..1: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
..1: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_row)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
object 19 7601.0 7677.2 -3781.5 7563.0
..1 20 7600.9 7681.0 -3780.5 7560.9 2.1304 1 0.1444
anova(fitfull,fitfull3) #There is strong evidence to use collumn however.
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
object: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
object: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench)
..1: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
..1: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
object 19 7601.0 7677.2 -3781.5 7563.0
..1 20 7597.9 7678.0 -3779.0 7557.9 5.13 1 0.02352 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#Modelling fixed effects.
fitfull3<-lmer(GM_TotalLeaf_Area~treatment*gluc_Conc*flav_Conc+BlackPathDam+WhiteFungDam+ThripsDam+Fern+(1|Family)+(1|gh_bench/gh_col), data=dat2)
#Removing three way interaction
fit.1<-update(fitfull3, ~.-treatment:gluc_Conc:flav_Conc)
anova(fitfull3,fit.1) #Good to remove
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
..1: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
..1: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
..1: treatment:gluc_Conc + treatment:flav_Conc + gluc_Conc:flav_Conc
object: GM_TotalLeaf_Area ~ treatment * gluc_Conc * flav_Conc + BlackPathDam +
object: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
..1 18 7600.9 7673.1 -3782.5 7564.9
object 20 7597.9 7678.0 -3779.0 7557.9 7.0377 2 0.02963 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
fit.2<-update(fit.1,~.-gluc_Conc:flav_Conc)
anova(fit.2,fit.1) #Good to remove.
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
fit.2: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.2: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
fit.2: treatment:gluc_Conc + treatment:flav_Conc
fit.1: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.1: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
fit.1: treatment:gluc_Conc + treatment:flav_Conc + gluc_Conc:flav_Conc
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.2 17 7601.7 7669.8 -3783.8 7567.7
fit.1 18 7600.9 7673.1 -3782.5 7564.9 2.7524 1 0.09711 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
fit.3<-update(fit.2,~.-treatment:flav_Conc)
anova(fit.2,fit.3) #Good to remove.
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
fit.3: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.3: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
fit.3: treatment:gluc_Conc
fit.2: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.2: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
fit.2: treatment:gluc_Conc + treatment:flav_Conc
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.3 15 7598.2 7658.3 -3784.1 7568.2
fit.2 17 7601.7 7669.8 -3783.8 7567.7 0.5123 2 0.774
fit.4<-update(fit.3,~.-treatment:gluc_Conc)
anova(fit.4,fit.3) #That significantly affected the predictive power. Did Not Remove
refitting model(s) with ML (instead of REML)
Data: dat2
Models:
fit.4: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.4: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col)
fit.3: GM_TotalLeaf_Area ~ treatment + gluc_Conc + flav_Conc + BlackPathDam +
fit.3: WhiteFungDam + ThripsDam + Fern + (1 | Family) + (1 | gh_bench/gh_col) +
fit.3: treatment:gluc_Conc
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.4 13 7602.4 7654.5 -3788.2 7576.4
fit.3 15 7598.2 7658.3 -3784.1 7568.2 8.1828 2 0.01672 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
fit.4<-update(fit.3,~.-flav_Conc)
#It seems as though flav conc has less sample size and so an anova cannot be done.I will refit the model with the same data set but without flav_Conc to test if it is significant.
fit.4<-lmer(GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
ThripsDam + Fern + treatment:gluc_Conc+ (1 | Family) + (1 | gh_bench/gh_col) ,data=dat2[!is.na(dat2$flav_Conc),])
fit.3<-lmer(GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
ThripsDam + Fern + treatment:gluc_Conc+flav_Conc+ (1 | Family) + (1 | gh_bench/gh_col) ,data=dat2[!is.na(dat2$flav_Conc),])
anova(fit.3,fit.4) #did not remove.
refitting model(s) with ML (instead of REML)
Data: dat2[!is.na(dat2$flav_Conc), ]
Models:
..1: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
..1: ThripsDam + Fern + treatment:gluc_Conc + (1 | Family) + (1 |
..1: gh_bench/gh_col)
object: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
object: ThripsDam + Fern + treatment:gluc_Conc + flav_Conc + (1 |
object: Family) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
..1 14 7600.5 7656.6 -3786.3 7572.5
object 15 7598.2 7658.3 -3784.1 7568.2 4.3304 1 0.03744 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#Flav_Conc is definitely a significant predictor on its own.
fit.4<-update(fit.3,~.-BlackPathDam)
anova(fit.3,fit.4) #do not remove black path dam.
refitting model(s) with ML (instead of REML)
Data: dat2[!is.na(dat2$flav_Conc), ]
Models:
..1: GM_TotalLeaf_Area ~ treatment + gluc_Conc + WhiteFungDam + ThripsDam +
..1: Fern + flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) +
..1: treatment:gluc_Conc
object: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
object: ThripsDam + Fern + treatment:gluc_Conc + flav_Conc + (1 |
object: Family) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
..1 14 7603.6 7659.7 -3787.8 7575.6
object 15 7598.2 7658.3 -3784.1 7568.2 7.4007 1 0.00652 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
fit.4<-update(fit.3,~.-WhiteFungDam)
anova(fit.3,fit.4) #removing whitefungdam.
refitting model(s) with ML (instead of REML)
Data: dat2[!is.na(dat2$flav_Conc), ]
Models:
..1: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + ThripsDam +
..1: Fern + flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) +
..1: treatment:gluc_Conc
object: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + WhiteFungDam +
object: ThripsDam + Fern + treatment:gluc_Conc + flav_Conc + (1 |
object: Family) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
..1 14 7596.4 7652.5 -3784.2 7568.4
object 15 7598.2 7658.3 -3784.1 7568.2 0.1653 1 0.6844
fit.5<-update(fit.4,~.-ThripsDam)
anova(fit.5,fit.4) #removing ThripsDam.
refitting model(s) with ML (instead of REML)
Data: dat2[!is.na(dat2$flav_Conc), ]
Models:
fit.5: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
fit.5: flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
fit.4: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + ThripsDam +
fit.4: Fern + flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) +
fit.4: treatment:gluc_Conc
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.5 13 7595.6 7647.7 -3784.8 7569.6
fit.4 14 7596.4 7652.5 -3784.2 7568.4 1.2195 1 0.2695
fit.6<-update(fit.5,~.-Fern)
anova(fit.6,fit.5) #Cannot Remove Fern.
refitting model(s) with ML (instead of REML)
Data: dat2[!is.na(dat2$flav_Conc), ]
Models:
fit.6: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + flav_Conc +
fit.6: (1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
fit.5: GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
fit.5: flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.6 12 7600.7 7648.8 -3788.4 7576.7
fit.5 13 7595.6 7647.7 -3784.8 7569.6 7.134 1 0.007564 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#fit.5 is therefore the best model
#testing to see if the effects remail using the whole data set, because flav_Conc removed 20 observations due to missing values.
fit.5.Whole<-lmer(GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
(1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc,data=dat2)
summary(fit.5.Whole)
Linear mixed model fit by REML t-tests use Satterthwaite
approximations to degrees of freedom [lmerMod]
Formula:
GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
(1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
Data: dat2
REML criterion at convergence: 7834.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.8150 -0.6034 0.0640 0.5573 3.5020
Random effects:
Groups Name Variance Std.Dev.
gh_col:gh_bench (Intercept) 358959 599.1
Family (Intercept) 102188 319.7
gh_bench (Intercept) 1905371 1380.4
Residual 6951563 2636.6
Number of obs: 426, groups:
gh_col:gh_bench, 28; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 13830.62 2107.60 110.96 6.562 1.76e-09 ***
treatmentgm -8384.71 2498.38 125.28 -3.356 0.00105 **
treatmentm -9186.17 2727.11 127.20 -3.368 0.00100 **
gluc_Conc -4026.29 1991.97 126.40 -2.021 0.04537 *
BlackPathDam -80.42 27.65 135.22 -2.909 0.00424 **
Fern -166.07 53.22 135.94 -3.120 0.00221 **
treatmentgm:gluc_Conc 4603.32 2514.14 125.19 1.831 0.06948 .
treatmentm:gluc_Conc 8356.75 2804.11 127.93 2.980 0.00345 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm glc_Cn BlckPD Fern trtmntg:_C
treatmentgm -0.761
treatmentm -0.698 0.592
gluc_Conc -0.944 0.795 0.729
BlackPathDm -0.078 0.036 0.083 0.036
Fern -0.008 0.010 -0.032 0.001 -0.143
trtmntgm:_C 0.745 -0.991 -0.581 -0.787 -0.047 -0.008
trtmntm:g_C 0.669 -0.567 -0.992 -0.707 -0.083 0.023 0.564
summary(fit.5)
Linear mixed model fit by REML ['lmerMod']
Formula:
GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
flav_Conc + (1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
Data: dat2[!is.na(dat2$flav_Conc), ]
REML criterion at convergence: 7443.5
Scaled residuals:
Min 1Q Median 3Q Max
-2.8937 -0.6083 0.0186 0.5802 3.4737
Random effects:
Groups Name Variance Std.Dev.
gh_col:gh_bench (Intercept) 457532 676.4
Family (Intercept) 151352 389.0
gh_bench (Intercept) 1901212 1378.8
Residual 6846335 2616.6
Number of obs: 406, groups:
gh_col:gh_bench, 28; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error t value
(Intercept) 12500.50 2173.11 5.752
treatmentgm -7461.78 2531.34 -2.948
treatmentm -9027.22 2816.04 -3.206
gluc_Conc -4727.86 2114.45 -2.236
BlackPathDam -73.25 28.03 -2.613
Fern -143.67 53.77 -2.672
flav_Conc 2289.97 1118.24 2.048
treatmentgm:gluc_Conc 3708.32 2548.62 1.455
treatmentm:gluc_Conc 8252.48 2909.58 2.836
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm glc_Cn BlckPD Fern flv_Cn
treatmentgm -0.767
treatmentm -0.688 0.590
gluc_Conc -0.841 0.741 0.668
BlackPathDm -0.092 0.041 0.093 -0.002
Fern -0.016 0.010 -0.049 -0.021 -0.140
flav_Conc -0.181 0.085 0.069 -0.281 0.113 0.065
trtmntgm:_C 0.749 -0.991 -0.579 -0.741 -0.050 -0.007 -0.066
trtmntm:g_C 0.659 -0.565 -0.992 -0.646 -0.093 0.042 -0.071
trtmntg:_C
treatmentgm
treatmentm
gluc_Conc
BlackPathDm
Fern
flav_Conc
trtmntgm:_C
trtmntm:g_C 0.562
#the results did not change. fit.5 is the best model.
Perhaps i need to include treatment in the nesting because each family may have different amount of individuals in different treatments, which could bias the fitness estimate of that family (the total estimate).
#Diagnostics #it didnt change the result at all and i think it makes more sense so i will keep it.
fit.5.Whole.tr<-lmer(GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
(1 | Family:treatment) + (1 | gh_bench/gh_col) + treatment:gluc_Conc,data=dat2)
plot(fit.5.Whole.tr)
summary(fit.5.Whole.tr)
Linear mixed model fit by REML t-tests use Satterthwaite
approximations to degrees of freedom [lmerMod]
Formula:
GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
(1 | Family:treatment) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
Data: dat2
REML criterion at convergence: 7835
Scaled residuals:
Min 1Q Median 3Q Max
-2.9387 -0.6079 0.0456 0.5633 3.4591
Random effects:
Groups Name Variance Std.Dev.
Family:treatment (Intercept) 162669 403.3
gh_col:gh_bench (Intercept) 378481 615.2
gh_bench (Intercept) 1893989 1376.2
Residual 6886655 2624.2
Number of obs: 426, groups:
Family:treatment, 69; gh_col:gh_bench, 28; gh_bench, 5
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 13801.70 2107.76 228.10 6.548 3.81e-10 ***
treatmentgm -8413.20 2506.43 448.30 -3.357 0.000856 ***
treatmentm -8921.50 2745.11 427.00 -3.250 0.001246 **
gluc_Conc -4009.11 1992.29 452.90 -2.012 0.044779 *
BlackPathDam -77.64 27.64 459.70 -2.809 0.005175 **
Fern -166.69 53.28 457.70 -3.129 0.001867 **
treatmentgm:gluc_Conc 4630.86 2519.57 451.70 1.838 0.066725 .
treatmentm:gluc_Conc 8071.45 2820.97 431.40 2.861 0.004425 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm glc_Cn BlckPD Fern trtmntg:_C
treatmentgm -0.762
treatmentm -0.699 0.586
gluc_Conc -0.944 0.794 0.728
BlackPathDm -0.076 0.034 0.079 0.034
Fern -0.009 0.011 -0.033 0.002 -0.140
trtmntgm:_C 0.745 -0.990 -0.575 -0.787 -0.044 -0.009
trtmntm:g_C 0.669 -0.561 -0.991 -0.707 -0.080 0.024 0.558
confint(fit.5.Whole)
Computing profile confidence intervals ...
2.5 % 97.5 %
.sig01 17.51414 1077.63029
.sig02 0.00000 722.19382
.sig03 582.90291 2871.46257
.sigma 2438.39069 2812.33052
(Intercept) 9768.13584 17915.80109
treatmentgm -13242.76168 -3507.00853
treatmentm -14520.80845 -3816.89927
gluc_Conc -7902.22126 -132.40434
BlackPathDam -135.12067 -26.45772
Fern -269.88524 -61.53812
treatmentgm:gluc_Conc -307.47814 9488.63612
treatmentm:gluc_Conc 2827.90704 13857.24758
coef(fit.5.Whole)
$`gh_col:gh_bench`
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam Fern
1:1 14457.00 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
1:2 13883.94 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
1:3 13885.42 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
1:4 13461.37 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2:1 13704.96 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2:2 13916.79 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2:3 14103.85 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2:4 14354.78 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2:5 13803.48 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
3:1 14208.50 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
3:2 14156.77 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
3:3 14323.58 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
3:4 14037.65 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4:1 13865.99 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4:2 14204.20 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4:3 14204.21 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4:4 13790.07 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4:5 13712.26 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5:1 13441.68 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5:2 13698.81 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5:3 13076.56 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5:4 13747.51 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5:5 13841.69 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
6:1 13688.24 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
6:2 13061.71 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
6:3 13459.58 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
6:4 13425.92 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
6:5 13740.94 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
treatmentgm:gluc_Conc treatmentm:gluc_Conc
1:1 4603.321 8356.747
1:2 4603.321 8356.747
1:3 4603.321 8356.747
1:4 4603.321 8356.747
2:1 4603.321 8356.747
2:2 4603.321 8356.747
2:3 4603.321 8356.747
2:4 4603.321 8356.747
2:5 4603.321 8356.747
3:1 4603.321 8356.747
3:2 4603.321 8356.747
3:3 4603.321 8356.747
3:4 4603.321 8356.747
4:1 4603.321 8356.747
4:2 4603.321 8356.747
4:3 4603.321 8356.747
4:4 4603.321 8356.747
4:5 4603.321 8356.747
5:1 4603.321 8356.747
5:2 4603.321 8356.747
5:3 4603.321 8356.747
5:4 4603.321 8356.747
5:5 4603.321 8356.747
6:1 4603.321 8356.747
6:2 4603.321 8356.747
6:3 4603.321 8356.747
6:4 4603.321 8356.747
6:5 4603.321 8356.747
$Family
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam
BWPEM1 13711.54 -8384.714 -9186.173 -4026.294 -80.42112
CBMCK1 13955.53 -8384.714 -9186.173 -4026.294 -80.42112
CRSOSO 14016.32 -8384.714 -9186.173 -4026.294 -80.42112
CWRIC2 13850.12 -8384.714 -9186.173 -4026.294 -80.42112
DVGM 13928.01 -8384.714 -9186.173 -4026.294 -80.42112
EFCC2 13738.55 -8384.714 -9186.173 -4026.294 -80.42112
JARI1 13930.65 -8384.714 -9186.173 -4026.294 -80.42112
JBBLB2 13613.39 -8384.714 -9186.173 -4026.294 -80.42112
JBCHY1 13478.37 -8384.714 -9186.173 -4026.294 -80.42112
JWBOY 14049.13 -8384.714 -9186.173 -4026.294 -80.42112
KVEDG1 13573.12 -8384.714 -9186.173 -4026.294 -80.42112
MAVBEL2 13854.25 -8384.714 -9186.173 -4026.294 -80.42112
MHBUR1 13879.65 -8384.714 -9186.173 -4026.294 -80.42112
MHNAT1 13823.37 -8384.714 -9186.173 -4026.294 -80.42112
MSMID1 13901.99 -8384.714 -9186.173 -4026.294 -80.42112
PDVRT1 13941.31 -8384.714 -9186.173 -4026.294 -80.42112
RULEB1 13871.36 -8384.714 -9186.173 -4026.294 -80.42112
SMAKC1 13826.50 -8384.714 -9186.173 -4026.294 -80.42112
SMITH1 14028.71 -8384.714 -9186.173 -4026.294 -80.42112
VRCAN 13767.34 -8384.714 -9186.173 -4026.294 -80.42112
VRPET2 13836.30 -8384.714 -9186.173 -4026.294 -80.42112
VSGARO1 13775.69 -8384.714 -9186.173 -4026.294 -80.42112
WSSWM3 13753.13 -8384.714 -9186.173 -4026.294 -80.42112
Fern treatmentgm:gluc_Conc treatmentm:gluc_Conc
BWPEM1 -166.0728 4603.321 8356.747
CBMCK1 -166.0728 4603.321 8356.747
CRSOSO -166.0728 4603.321 8356.747
CWRIC2 -166.0728 4603.321 8356.747
DVGM -166.0728 4603.321 8356.747
EFCC2 -166.0728 4603.321 8356.747
JARI1 -166.0728 4603.321 8356.747
JBBLB2 -166.0728 4603.321 8356.747
JBCHY1 -166.0728 4603.321 8356.747
JWBOY -166.0728 4603.321 8356.747
KVEDG1 -166.0728 4603.321 8356.747
MAVBEL2 -166.0728 4603.321 8356.747
MHBUR1 -166.0728 4603.321 8356.747
MHNAT1 -166.0728 4603.321 8356.747
MSMID1 -166.0728 4603.321 8356.747
PDVRT1 -166.0728 4603.321 8356.747
RULEB1 -166.0728 4603.321 8356.747
SMAKC1 -166.0728 4603.321 8356.747
SMITH1 -166.0728 4603.321 8356.747
VRCAN -166.0728 4603.321 8356.747
VRPET2 -166.0728 4603.321 8356.747
VSGARO1 -166.0728 4603.321 8356.747
WSSWM3 -166.0728 4603.321 8356.747
$gh_bench
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam Fern
1 15861.65 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
2 13504.05 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
3 14199.31 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
4 12947.15 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
5 12640.96 -8384.714 -9186.173 -4026.294 -80.42112 -166.0728
treatmentgm:gluc_Conc treatmentm:gluc_Conc
1 4603.321 8356.747
2 4603.321 8356.747
3 4603.321 8356.747
4 4603.321 8356.747
5 4603.321 8356.747
attr(,"class")
[1] "coef.mer"
#no heteroscedasticity
plot(fit.5.Whole)
#fairly normal.
qqnorm(resid(fit.5.Whole))
#Flavonoid is significant, but just barely, there are also alot of samples missing there and I know there is a decent amount of sampling variance there, so i think i will ignore the small positive eeffect of flavonoids on the reduced dataset and continue with the whole one.
#Conclusion, there is a benefit to producing glucosniolates in the maple treamtment and it is detrimental in the alone treatment. In the garlic mustard treatment, there were mixed results.
#Visualization
library(ggplot2)
coef(fit.5.Whole.tr)
$`Family:treatment`
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam
BWPEM1:a 13895.17 -8413.2 -8921.5 -4009.115 -77.64257
BWPEM1:gm 13618.16 -8413.2 -8921.5 -4009.115 -77.64257
BWPEM1:m 13690.67 -8413.2 -8921.5 -4009.115 -77.64257
CBMCK1:a 13788.42 -8413.2 -8921.5 -4009.115 -77.64257
CBMCK1:gm 13955.57 -8413.2 -8921.5 -4009.115 -77.64257
CBMCK1:m 13878.20 -8413.2 -8921.5 -4009.115 -77.64257
CRSOSO:a 13733.00 -8413.2 -8921.5 -4009.115 -77.64257
CRSOSO:gm 13882.56 -8413.2 -8921.5 -4009.115 -77.64257
CRSOSO:m 14130.23 -8413.2 -8921.5 -4009.115 -77.64257
CWRIC2:a 13779.46 -8413.2 -8921.5 -4009.115 -77.64257
CWRIC2:gm 13735.72 -8413.2 -8921.5 -4009.115 -77.64257
CWRIC2:m 13925.59 -8413.2 -8921.5 -4009.115 -77.64257
DVGM:a 13684.40 -8413.2 -8921.5 -4009.115 -77.64257
DVGM:gm 13859.63 -8413.2 -8921.5 -4009.115 -77.64257
DVGM:m 14039.02 -8413.2 -8921.5 -4009.115 -77.64257
EFCC2:a 13707.57 -8413.2 -8921.5 -4009.115 -77.64257
EFCC2:gm 13740.23 -8413.2 -8921.5 -4009.115 -77.64257
EFCC2:m 13794.89 -8413.2 -8921.5 -4009.115 -77.64257
JARI1:a 13907.50 -8413.2 -8921.5 -4009.115 -77.64257
JARI1:gm 13850.78 -8413.2 -8921.5 -4009.115 -77.64257
JARI1:m 13826.83 -8413.2 -8921.5 -4009.115 -77.64257
JBBLB2:a 13827.29 -8413.2 -8921.5 -4009.115 -77.64257
JBBLB2:gm 13763.45 -8413.2 -8921.5 -4009.115 -77.64257
JBBLB2:m 13418.74 -8413.2 -8921.5 -4009.115 -77.64257
JBCHY1:a 13773.18 -8413.2 -8921.5 -4009.115 -77.64257
JBCHY1:gm 13547.51 -8413.2 -8921.5 -4009.115 -77.64257
JBCHY1:m 13444.65 -8413.2 -8921.5 -4009.115 -77.64257
JWBOY:a 13881.99 -8413.2 -8921.5 -4009.115 -77.64257
JWBOY:gm 14127.63 -8413.2 -8921.5 -4009.115 -77.64257
JWBOY:m 13783.77 -8413.2 -8921.5 -4009.115 -77.64257
KVEDG1:a 13591.14 -8413.2 -8921.5 -4009.115 -77.64257
KVEDG1:gm 13670.93 -8413.2 -8921.5 -4009.115 -77.64257
KVEDG1:m 13673.42 -8413.2 -8921.5 -4009.115 -77.64257
MAVBEL2:a 13807.00 -8413.2 -8921.5 -4009.115 -77.64257
MAVBEL2:gm 13801.01 -8413.2 -8921.5 -4009.115 -77.64257
MAVBEL2:m 13846.82 -8413.2 -8921.5 -4009.115 -77.64257
MHBUR1:a 13989.12 -8413.2 -8921.5 -4009.115 -77.64257
MHBUR1:gm 13753.58 -8413.2 -8921.5 -4009.115 -77.64257
MHBUR1:m 13748.85 -8413.2 -8921.5 -4009.115 -77.64257
MHNAT1:a 13812.01 -8413.2 -8921.5 -4009.115 -77.64257
MHNAT1:gm 13744.95 -8413.2 -8921.5 -4009.115 -77.64257
MHNAT1:m 13835.96 -8413.2 -8921.5 -4009.115 -77.64257
MSMID1:a 13896.21 -8413.2 -8921.5 -4009.115 -77.64257
MSMID1:gm 13829.53 -8413.2 -8921.5 -4009.115 -77.64257
MSMID1:m 13805.22 -8413.2 -8921.5 -4009.115 -77.64257
PDVRT1:a 13730.17 -8413.2 -8921.5 -4009.115 -77.64257
PDVRT1:gm 13945.93 -8413.2 -8921.5 -4009.115 -77.64257
PDVRT1:m 13922.60 -8413.2 -8921.5 -4009.115 -77.64257
RULEB1:a 13873.96 -8413.2 -8921.5 -4009.115 -77.64257
RULEB1:gm 13751.97 -8413.2 -8921.5 -4009.115 -77.64257
RULEB1:m 13850.02 -8413.2 -8921.5 -4009.115 -77.64257
SMAKC1:a 13726.88 -8413.2 -8921.5 -4009.115 -77.64257
SMAKC1:gm 13725.68 -8413.2 -8921.5 -4009.115 -77.64257
SMAKC1:m 13946.78 -8413.2 -8921.5 -4009.115 -77.64257
SMITH1:a 13669.61 -8413.2 -8921.5 -4009.115 -77.64257
SMITH1:gm 14089.80 -8413.2 -8921.5 -4009.115 -77.64257
SMITH1:m 13997.47 -8413.2 -8921.5 -4009.115 -77.64257
VRCAN:a 13813.62 -8413.2 -8921.5 -4009.115 -77.64257
VRCAN:gm 13783.00 -8413.2 -8921.5 -4009.115 -77.64257
VRCAN:m 13692.68 -8413.2 -8921.5 -4009.115 -77.64257
VRPET2:a 13874.44 -8413.2 -8921.5 -4009.115 -77.64257
VRPET2:gm 13859.00 -8413.2 -8921.5 -4009.115 -77.64257
VRPET2:m 13685.91 -8413.2 -8921.5 -4009.115 -77.64257
VSGARO1:a 13841.03 -8413.2 -8921.5 -4009.115 -77.64257
VSGARO1:gm 13857.35 -8413.2 -8921.5 -4009.115 -77.64257
VSGARO1:m 13612.91 -8413.2 -8921.5 -4009.115 -77.64257
WSSWM3:a 13836.05 -8413.2 -8921.5 -4009.115 -77.64257
WSSWM3:gm 13545.27 -8413.2 -8921.5 -4009.115 -77.64257
WSSWM3:m 13887.97 -8413.2 -8921.5 -4009.115 -77.64257
Fern treatmentgm:gluc_Conc treatmentm:gluc_Conc
BWPEM1:a -166.6875 4630.863 8071.455
BWPEM1:gm -166.6875 4630.863 8071.455
BWPEM1:m -166.6875 4630.863 8071.455
CBMCK1:a -166.6875 4630.863 8071.455
CBMCK1:gm -166.6875 4630.863 8071.455
CBMCK1:m -166.6875 4630.863 8071.455
CRSOSO:a -166.6875 4630.863 8071.455
CRSOSO:gm -166.6875 4630.863 8071.455
CRSOSO:m -166.6875 4630.863 8071.455
CWRIC2:a -166.6875 4630.863 8071.455
CWRIC2:gm -166.6875 4630.863 8071.455
CWRIC2:m -166.6875 4630.863 8071.455
DVGM:a -166.6875 4630.863 8071.455
DVGM:gm -166.6875 4630.863 8071.455
DVGM:m -166.6875 4630.863 8071.455
EFCC2:a -166.6875 4630.863 8071.455
EFCC2:gm -166.6875 4630.863 8071.455
EFCC2:m -166.6875 4630.863 8071.455
JARI1:a -166.6875 4630.863 8071.455
JARI1:gm -166.6875 4630.863 8071.455
JARI1:m -166.6875 4630.863 8071.455
JBBLB2:a -166.6875 4630.863 8071.455
JBBLB2:gm -166.6875 4630.863 8071.455
JBBLB2:m -166.6875 4630.863 8071.455
JBCHY1:a -166.6875 4630.863 8071.455
JBCHY1:gm -166.6875 4630.863 8071.455
JBCHY1:m -166.6875 4630.863 8071.455
JWBOY:a -166.6875 4630.863 8071.455
JWBOY:gm -166.6875 4630.863 8071.455
JWBOY:m -166.6875 4630.863 8071.455
KVEDG1:a -166.6875 4630.863 8071.455
KVEDG1:gm -166.6875 4630.863 8071.455
KVEDG1:m -166.6875 4630.863 8071.455
MAVBEL2:a -166.6875 4630.863 8071.455
MAVBEL2:gm -166.6875 4630.863 8071.455
MAVBEL2:m -166.6875 4630.863 8071.455
MHBUR1:a -166.6875 4630.863 8071.455
MHBUR1:gm -166.6875 4630.863 8071.455
MHBUR1:m -166.6875 4630.863 8071.455
MHNAT1:a -166.6875 4630.863 8071.455
MHNAT1:gm -166.6875 4630.863 8071.455
MHNAT1:m -166.6875 4630.863 8071.455
MSMID1:a -166.6875 4630.863 8071.455
MSMID1:gm -166.6875 4630.863 8071.455
MSMID1:m -166.6875 4630.863 8071.455
PDVRT1:a -166.6875 4630.863 8071.455
PDVRT1:gm -166.6875 4630.863 8071.455
PDVRT1:m -166.6875 4630.863 8071.455
RULEB1:a -166.6875 4630.863 8071.455
RULEB1:gm -166.6875 4630.863 8071.455
RULEB1:m -166.6875 4630.863 8071.455
SMAKC1:a -166.6875 4630.863 8071.455
SMAKC1:gm -166.6875 4630.863 8071.455
SMAKC1:m -166.6875 4630.863 8071.455
SMITH1:a -166.6875 4630.863 8071.455
SMITH1:gm -166.6875 4630.863 8071.455
SMITH1:m -166.6875 4630.863 8071.455
VRCAN:a -166.6875 4630.863 8071.455
VRCAN:gm -166.6875 4630.863 8071.455
VRCAN:m -166.6875 4630.863 8071.455
VRPET2:a -166.6875 4630.863 8071.455
VRPET2:gm -166.6875 4630.863 8071.455
VRPET2:m -166.6875 4630.863 8071.455
VSGARO1:a -166.6875 4630.863 8071.455
VSGARO1:gm -166.6875 4630.863 8071.455
VSGARO1:m -166.6875 4630.863 8071.455
WSSWM3:a -166.6875 4630.863 8071.455
WSSWM3:gm -166.6875 4630.863 8071.455
WSSWM3:m -166.6875 4630.863 8071.455
$`gh_col:gh_bench`
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam Fern
1:1 14467.92 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
1:2 13849.13 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
1:3 13862.44 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
1:4 13439.49 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2:1 13683.40 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2:2 13879.28 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2:3 14072.10 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2:4 14358.82 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2:5 13763.75 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
3:1 14201.94 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
3:2 14162.16 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
3:3 14323.11 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
3:4 14002.27 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4:1 13844.55 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4:2 14190.87 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4:3 14181.11 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4:4 13753.79 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4:5 13678.42 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5:1 13379.44 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5:2 13649.49 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5:3 13031.63 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5:4 13701.36 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5:5 13811.83 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
6:1 13639.40 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
6:2 13012.93 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
6:3 13409.28 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
6:4 13379.28 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
6:5 13718.56 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
treatmentgm:gluc_Conc treatmentm:gluc_Conc
1:1 4630.863 8071.455
1:2 4630.863 8071.455
1:3 4630.863 8071.455
1:4 4630.863 8071.455
2:1 4630.863 8071.455
2:2 4630.863 8071.455
2:3 4630.863 8071.455
2:4 4630.863 8071.455
2:5 4630.863 8071.455
3:1 4630.863 8071.455
3:2 4630.863 8071.455
3:3 4630.863 8071.455
3:4 4630.863 8071.455
4:1 4630.863 8071.455
4:2 4630.863 8071.455
4:3 4630.863 8071.455
4:4 4630.863 8071.455
4:5 4630.863 8071.455
5:1 4630.863 8071.455
5:2 4630.863 8071.455
5:3 4630.863 8071.455
5:4 4630.863 8071.455
5:5 4630.863 8071.455
6:1 4630.863 8071.455
6:2 4630.863 8071.455
6:3 4630.863 8071.455
6:4 4630.863 8071.455
6:5 4630.863 8071.455
$gh_bench
(Intercept) treatmentgm treatmentm gluc_Conc BlackPathDam Fern
1 15835.50 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
2 13469.57 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
3 14149.24 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
4 12924.83 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
5 12629.37 -8413.2 -8921.5 -4009.115 -77.64257 -166.6875
treatmentgm:gluc_Conc treatmentm:gluc_Conc
1 4630.863 8071.455
2 4630.863 8071.455
3 4630.863 8071.455
4 4630.863 8071.455
5 4630.863 8071.455
attr(,"class")
[1] "coef.mer"
summary(fit.5.Whole)
Linear mixed model fit by REML t-tests use Satterthwaite
approximations to degrees of freedom [lmerMod]
Formula:
GM_TotalLeaf_Area ~ treatment + gluc_Conc + BlackPathDam + Fern +
(1 | Family) + (1 | gh_bench/gh_col) + treatment:gluc_Conc
Data: dat2
REML criterion at convergence: 7834.8
Scaled residuals:
Min 1Q Median 3Q Max
-2.8150 -0.6034 0.0640 0.5573 3.5020
Random effects:
Groups Name Variance Std.Dev.
gh_col:gh_bench (Intercept) 358959 599.1
Family (Intercept) 102188 319.7
gh_bench (Intercept) 1905371 1380.4
Residual 6951563 2636.6
Number of obs: 426, groups:
gh_col:gh_bench, 28; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 13830.62 2107.60 110.96 6.562 1.76e-09 ***
treatmentgm -8384.71 2498.38 125.28 -3.356 0.00105 **
treatmentm -9186.17 2727.11 127.20 -3.368 0.00100 **
gluc_Conc -4026.29 1991.97 126.40 -2.021 0.04537 *
BlackPathDam -80.42 27.65 135.22 -2.909 0.00424 **
Fern -166.07 53.22 135.94 -3.120 0.00221 **
treatmentgm:gluc_Conc 4603.32 2514.14 125.19 1.831 0.06948 .
treatmentm:gluc_Conc 8356.75 2804.11 127.93 2.980 0.00345 **
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm glc_Cn BlckPD Fern trtmntg:_C
treatmentgm -0.761
treatmentm -0.698 0.592
gluc_Conc -0.944 0.795 0.729
BlackPathDm -0.078 0.036 0.083 0.036
Fern -0.008 0.010 -0.032 0.001 -0.143
trtmntgm:_C 0.745 -0.991 -0.581 -0.787 -0.047 -0.008
trtmntm:g_C 0.669 -0.567 -0.992 -0.707 -0.083 0.023 0.564
ggplot(dat3)+
geom_point(aes(y=GM_TotalLeaf_Area,x=gluc_Conc,colour=treatment))+
geom_abline(intercept=13830.62,slope=-4026.29,colour="red")+
geom_abline(intercept=13830.62-9186.17,slope=-4026.29+8356.75,colour="blue")
ggplot(dat3)+
geom_point(aes(y=GM_TotalLeaf_Area,x=BlackPathDam,colour=treatment))
ggplot(dat2)+
geom_point(aes(y=GM_TotalLeaf_Area,x=Fern))+
geom_abline(intercept=13830.62,slope=-166.07)
We wouldnt really expect to have the power to detect a GXE effect on fitness here, because each genotype will have a specific degree of gluc concentration. I dont think that i can do a model using the variables from leaf data, because the dependent variable is the same for each leaf, they come from the same plant.
#——————————————————————- GROWTHRATE #——————————————————————-
#Repeat this analysis on RGR1
GrowthDatMean<-left_join(dat3,GrowthDatMean,by=c("Family","treatment"))
Column `Family` joining character vector and factor, coercing into character vectorColumn `treatment` joining factors with different levels, coercing to character vector
#Modelling relative growth rate 1. Pathogens and fern abundance was removed from this analysis because they did not appear this early in the experiment.
#Modelling fixed effects. Ignoring flavonoids for the time being.
fitfull3<-lmer(RGR1~treatment*gluc_Conc+(1|Family:treatment)+(1|gh_bench/gh_col), data=GrowthDat)
#Removing three way interaction
fit.1<-update(fitfull3, ~.-gluc_Conc:treatment)
anova(fitfull3,fit.1) #Good to remove
fit.2<-update(fit.1,~.-gluc_Conc)
anova(fit.2,fit.1) #Good to remove.
fit.3<-update(fit.2,~.-treatment)
anova(fit.2,fit.3) #Significant
Treament was the only predictor for growth rate 1.
#Modelling on growth rate 2.
anova(fit.2,fit.3) #Significant
refitting model(s) with ML (instead of REML)
Data: GrowthDat
Models:
fit.3: RGR2 ~ (1 | Family:treatment) + (1 | gh_bench/gh_col)
fit.2: RGR2 ~ treatment + (1 | Family:treatment) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.3 5 -11315 -11295 5662.6 -11325
fit.2 7 -11331 -11303 5672.7 -11345 20.214 2 4.079e-05 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Treament was the only predictor for growth rate 2.
#Modelling on growth rate 3.
anova(fit.6,fit.7) #Not Significant
refitting model(s) with ML (instead of REML)
Data: GrowthDat[!is.na(GrowthDat$ThripsDam), ]
Models:
fit.7: RGR3 ~ (1 | Family:treatment) + (1 | gh_bench/gh_col)
fit.6: RGR3 ~ Fern + (1 | Family:treatment) + (1 | gh_bench/gh_col)
Df AIC BIC logLik deviance Chisq Chi Df Pr(>Chisq)
fit.7 5 -11574 -11554 5792.0 -11584
fit.6 6 -11575 -11551 5793.6 -11587 3.1064 1 0.07799 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#Modelling on growth rate 4. This one include the other variables as they were now involved in the experiement around this time.
summary(fit.2)
Linear mixed model fit by REML ['lmerMod']
Formula:
RGR4 ~ treatment * gluc_Conc + Fern + BlackPathDam + (1 | Family:treatment)
Data: GrowthDat
REML criterion at convergence: -11683
Scaled residuals:
Min 1Q Median 3Q Max
-4.0620 -0.4092 0.0438 0.5709 5.4749
Random effects:
Groups Name Variance Std.Dev.
Family:treatment (Intercept) 2.418e-16 1.555e-08
Residual 1.090e-14 1.044e-07
Number of obs: 408, groups: Family:treatment, 69
Fixed effects:
Estimate Std. Error t value
(Intercept) -1.595e-07 8.234e-08 -1.937
treatmentgm 2.462e-07 1.021e-07 2.410
treatmentm 2.626e-07 1.111e-07 2.363
gluc_Conc 1.715e-07 8.105e-08 2.116
Fern 4.398e-09 1.991e-09 2.209
BlackPathDam -5.202e-09 1.105e-09 -4.709
treatmentgm:gluc_Conc -2.345e-07 1.026e-07 -2.285
treatmentm:gluc_Conc -3.034e-07 1.138e-07 -2.665
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm glc_Cn Fern BlckPD trtmntg:_C
treatmentgm -0.804
treatmentm -0.744 0.597
gluc_Conc -0.992 0.799 0.736
Fern -0.005 0.000 -0.046 -0.005
BlackPathDm -0.078 0.035 0.100 0.034 -0.102
trtmntgm:_C 0.785 -0.991 -0.584 -0.790 0.006 -0.048
trtmntm:g_C 0.712 -0.571 -0.992 -0.715 0.037 -0.100 0.567
#The next step to support this analysis would be to determine if glucosinolate content explains competitive ability against maples.
#Look at if there is a negative correlation between biomass in the maple and alone treatment.
#—————————————– #Part 2
#assessing
#Significant testing (doing seperately because flav conc and gluc conc are correlated)
#WhitePathDam
#BlackPathDam
summary(fit.4)
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: poisson ( log )
Formula: BlackPathDam ~ flav_Conc + (1 | Family/Tag) + (1 | gh_bench)
Data: dat[!is.na(dat$flav_Conc), ]
AIC BIC logLik deviance df.resid
3504.2 3526.8 -1747.1 3494.2 672
Scaled residuals:
Min 1Q Median 3Q Max
-3.0421 -0.4065 -0.2289 0.0876 4.0742
Random effects:
Groups Name Variance Std.Dev.
Tag:Family (Intercept) 0.971171 0.98548
Family (Intercept) 0.052429 0.22897
gh_bench (Intercept) 0.003662 0.06051
Number of obs: 677, groups: Tag:Family, 423; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.6816 0.1579 10.651 < 2e-16 ***
flav_Conc -1.1150 0.1644 -6.781 1.19e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
flav_Conc -0.856
#ThripsDam
summary(fit.4) #Flav conc is again, a very significant predictor.
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: poisson ( log )
Formula: BlackPathDam ~ flav_Conc + (1 | Family/Tag) + (1 | gh_bench)
Data: dat[!is.na(dat$flav_Conc), ]
AIC BIC logLik deviance df.resid
3504.2 3526.8 -1747.1 3494.2 672
Scaled residuals:
Min 1Q Median 3Q Max
-3.0421 -0.4065 -0.2289 0.0876 4.0742
Random effects:
Groups Name Variance Std.Dev.
Tag:Family (Intercept) 0.971171 0.98548
Family (Intercept) 0.052429 0.22897
gh_bench (Intercept) 0.003662 0.06051
Number of obs: 677, groups: Tag:Family, 423; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 1.6816 0.1579 10.651 < 2e-16 ***
flav_Conc -1.1150 0.1644 -6.781 1.19e-11 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr)
flav_Conc -0.856
summary(fit.1) #Flav conc and treatment predict fern abudnance.
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: poisson ( log )
Formula: Fern ~ treatment + flav_Conc + (1 | Family) + (1 | gh_bench/gh_col)
Data: dat[!is.na(dat$flav_Conc), ]
AIC BIC logLik deviance df.resid
1804.4 1836.1 -895.2 1790.4 680
Scaled residuals:
Min 1Q Median 3Q Max
-2.5073 -0.6085 -0.2213 -0.1154 16.0269
Random effects:
Groups Name Variance Std.Dev.
gh_col:gh_bench (Intercept) 3.636 1.907
Family (Intercept) 1.701 1.304
gh_bench (Intercept) 1.153 1.074
Number of obs: 687, groups: gh_col:gh_bench, 28; Family, 23; gh_bench, 5
Fixed effects:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -2.1683 0.7522 -2.883 0.003944 **
treatmentgm -0.2786 0.1498 -1.859 0.063030 .
treatmentm 0.4068 0.1122 3.625 0.000289 ***
flav_Conc -0.6571 0.2737 -2.401 0.016365 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation of Fixed Effects:
(Intr) trtmntg trtmntm
treatmentgm -0.148
treatmentm -0.165 0.500
flav_Conc -0.325 0.221 0.252
#effect on maples